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Abstract 



A dilute Bose system with Bose-Einstein condensate is considered. It is 
shown that the Hartree-Fock-Bogolubov approximation can be made both 
conserving as well as gapless. This is achieved by taking into account all 
physical normalization conditions, that is, the normalization condition for the 
condensed particles and that for the total number of particles. Two Lagrange 
multipliers, introduced for preserving these normalization conditions, make 
the consideration completely self-consistent. 
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I. INTRODUCTION 



The physical properties of Bose gases exhibiting Bose-Einstein condensation are now a 
topic of intensive research, both experimental and theoretical [1-5]. The physics of dilute 
weakly interacting Bose gas has been studied especially well using Bogoliubov's model and 
operator techniques [6-9]. Despite the apparent simplicity of the system, some principal 
problems in its theoretical description have not yet been well understood. Most notorious is 
the so-called dilemma of conserving versus gapless approximations, for which we shall offer 
a possible solution in this paper. 

The above dilemma arises, when one attempts to describe the weakly interacting Bose 
gas at finite (not asymptotically small) diluteness parameter. It is absent in the original 
Bogoliubov approximation [6-9], which is valid only at very low temperatures T — > and 
asymptotically weak interaction. The first natural extension to finite temperatures and 
diluteness parameters is a bosonic version of the Hartree-Fock-Bogoliubov (HFB) approxi- 
mation. This is a self-consistent approximation which is guaranteed to respect all conser- 
vation laws that follow from the underlying symmetries of the Hamiltonian. However, this 
approximation turns out to render a gap in the spectrum of collective excitations in the 
condensed phase in which the global U(l) gauge symmetry of the theory is broken. Also the 
Girardeau- Arnowitt approximation [10] displays this gap, since it is equivalent to the HFB 
approximation. This contradicts the fundamental theorems of Hugenholtz and Pines [11], 
Goldstone [12], and Bogoliubov [9], according to which the spectrum of collective excitations 
in the symmetry-broken phase has to be gapless. It is well known [13-15] that the HFB 
approximation can be formulated as a variational approximation. There exist several other 
variational approximations, the so-called ^-derivable approximations or higher effective ac- 
tions [16-20], which also respect the conservation laws but lead to a gap in the spectrum. 
When trying to remove the gap, one usually violates the conservation laws and runs into 
other thermodynamic inconsistencies. The various approximations that have been studied 
are typically classified as either conserving or gapless. This classification and the related 
dilemma were first emphasized by Hohenberg and Martin [21], and later discussed in many 
publications, for instance, by Baym and Grinstein [19]. A very detailed discussion of this 
problem, with many citations, has recently been done by Andersen [4]. 

In order to remove the gap in the HFB approximation, one often invokes a trick of 
neglecting the anomalous averages, calling this the "Popov approximation". However, a 
glance at the original works by Popov [22-25] shows that he has never suggested this trick. 
What he actually considered was a narrow region of temperatures T in the vicinity of the 
condensation temperature T c . When T — > T c , then the condensate density tends to zero. 
The anomalous averages, being proportional to the condensate density, tend to zero together 
with the latter, when T — > T c . As a result, their contribution becomes automatically small, 
without any special assumptions. Far below T c , however, the anomalous averages can be 
very large, and Popov has never proposed to neglect them there. It is straightforward to 
demonstrate by direct calculations that at low temperatures T < T c the anomalous averages 
become of the same order as the normal averages. They can even be much larger than the 
latter [26], so that Popov would certainly not have been keen on proposing the so-called 
" Popov approximation" . 

Moreover, preserving the anomalous averages for the phase with broken U(l) gauge 
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symmetry is principally important, since their negligence makes the system unstable [26,27]. 

There have been several attempts to make the HFB approximation gapless by adding 
to this approximation some higher-order terms involving the Bethe-Salpeter or T-matrix 
approximations. The idea to modify in this way the mean-field HFB approximation goes 
back to Kirzhnits and Linde [28] and Baym and Grinstein [19]. Several such modifications of 
the HFB approximation have been considered, in which additional terms are either motivated 
by higher-order approximations [29,30] or just added phenomenologically [31,32]. This type 
of modifications [28-32] possesses a number of deficiencies, which do not permit to accept 
this approach as a solution of the problem. A good analysis of these deficiencies has already 
been done by Baym and Grinstein [19] and recently by Andersen [4]. 

First of all, the method of adding to the theory some additional terms in order to cancel 
the spectrum gap is ambiguous, not following from a general physical law. As a result, the 
additional terms are not uniquely defined. There are no general grounds to decide which of 
the variants is better. 

Second, the way of mixing different approximations is not self-consistent. This is what 
Bogoliubov [9] called the "mismatch of approximations". Although such a mismatch can 
make the gap disappear, it is usually inconsistent with some thermodynamic equations. For 
instance, the chemical potential defined by the Hugenholtz-Pines theorem to yield a gapless 
spectrum, does not coincide with the chemical potential found from the minimization of 
the thermodynamic potential with respect to the condensate density. This discrepancy is a 
general feature of all non-self-consistent approximations, due to which they cannot properly 
be called gapless [4]. 

Moreover, if the price of making the spectrum gapless is that the thermodynamic po- 
tential cannot be minimized, this implies that the system becomes thermodynamically un- 
stable. Thus one has the unpleasant alternative: either the system is gapless but unstable 
or seemingly stable but with a gap. This is a particular case of the general problem of 
thermodynamic self-consistency [33,34]. 

Problems with the thermodynamics of the model have the unpleasant consequences of 
modifying the order of the phase transition [33,34]. This happens for all approaches with 
mismatched approximations. Baym and Grinstein [19] emphasized that attempts to modify 
the HFB approximation in a non-self-consistent way lead to a first-order condensation tran- 
sition, instead of the observed second-order one. Further works [4,31,32] confirmed that it 
is a common feature of all non-self-consistent mean-field approximations. 

In the present paper, we show that the HFB approximation can be made both conserving 
and gapless, while avoiding the mismatch of approximations, thus being completely self- 
consistent. The solution of the problem is possible by taking into account two existing 
normalization conditions instead of one in all previous approximations. This makes the 
HFB approximation self-consistent and gapless, without any tricks or additional terms. As 
a consequence, condensation transition becomes second order, as it should be. 

Throughout the paper, the system of units is employed with Planck and Boltzmann 
constants equal to unity, H = 1, ks = 1- 
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II. GRAND-CANONICAL HAMILTONIAN 



We consider a dilute Bose gas, whose particle interactions are modeled by the contact 
potential 

V(x) = g6(x) , (1) 

with the interaction intensity 

9 = ^- (2) 
m 

expressed through the s-wave scattering length a s and particle mass m. The energy operator 
reads 



(3) 



with the field operators satisfying the Bose commutation relations. The operator 

counting the total number of particles is 

N = J rf 3 x^ t (x)^(x). (4) 

In what follows, we treat a uniform system in thermodynamic equilibrium, since this is the 
simplest situation where the above-discussed problems reveal themselves. 

If the temperature in a Bose system falls below the condensation temperature T c , U(l) 
gauge symmetry becomes broken. The symmetry breaking is taken into account by the 
Bogoliubov shift [8,9] of the field operator 

^(x) — ► tt + ^i(x), (5) 

where \1/ is the condensate order parameter, which in uniform systems is independent of 
x, and is the field operator of the uncondensed particles, satisfying the same Bose 

commutation relations as ^(x). Another method of first separating the zero- momentum 
components of the field and then replacing them by classical numbers [6,7] is completely 
equivalent to the shift (5), as has been rigorously proved by Ginibre [35]. The field , 0i( x ) 
has no zero-momentum component so that 

0Mx)> = (6) 
and ^ and ipi (x) are orthogonal to each other: 

yV-0!(x) d 3 x = 0. (7) 
The condensate order parameter \l> defines the density of condensed particles 

Po = l^l 2 • (8) 

Above T c where ^ vanishes, there is no condensate. Below T c , one has ^ ^ and thus a 
finite condensate density p . 
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The free energy of the system is 

F = -ThxTr e~ pA , (9) 

with (3 = 1/T. After substituting shift (5) into H, the model contains two field variables, the 
condensate order parameter \1/ and the space-dependent field ^i(x). They give rise to two 
normalization conditions. One is related to the definition of the total number of particles 

N — (N) . (10) 

A second normalization condition fixes the number of condensed particles 

N = poV = \V\ 2 V . (11) 

In stable equilibrium, the free energy gains a minimum under the normalization con- 
ditions (10) and (11). This conditional minimization is equivalent to the unconditional 
minimum of the grand-canonical potential 

D, = -TlnTr e~ 0H , (12) 

with the grand-canonical Hamiltonian 

H = H-n Q N Q -nN , (13) 

in which /i and fi are the Lagrange multipliers enforcing the normalization conditions (10) 
and (11). The minimum of the grand-canonical potential (12) is determined by the equations 

» ^72 > • (14) 



dN dN$ 
The first derivative is given by the expectation value 

dfl I dH \ 
dN~ ~ \dN~o/ • 

The second derivative is calculated from 



(15) 



d 2 n ld 2 H\ nA2 fdH\ 

+ P A2 {^)> (16) 



dNg \dN^/ ' ^ \dN Q 
with the notation 

A 2 (6) = (O 2 ) - (6) 2 

for the dispersion of an operator O. 

For a uniform system, the field operator of the uncondensed particles is expandable in 
plane waves as 

iMx) = j>^. (17) 

k^o v V 



5 



Performing the Bogoliubov shift (5), together with the expansion (17), we obtain for the 
grand-canonical Hamiltonian (13) the sum of five terms 

4 

H = E # (n) , (18) 

n=0 

depending on the the power of the operators ^i(x). The zero-order term 

if(°)= (lp g-s^N , (19) 

with 

e = fio + fi, (20) 

is free of V'i(x). The first-order term in Eq. (19) vanishes, since the decomposition (17) 
contains only nonzero momenta. The second-order term is 

H (2) = E ( f- + 2 Po9 ~ 4«k + P ~f E (4« f - k + a- k a k ) . (21) 

In the third-order term 

^ ^ = \fv 9 ^ K^P+^-P + a - P «P +q «q) , (22) 

p,q 

the prime on the summation symbol implies that p ^ 0, q ^ 0, and p + q ^ 0. In the 
fourth-order term 

ff(4) 4EE'«iW- (23) 

ZV k P,q 

the summation does not include any zero-momentum operators, so that the prime means 
that p^0,q^0,k+p^0, and k - q ^ 0. 

The field operators in momentum space ak satisfy the following conditions. From Eq. 
(7) we have 

(a k > = (k ^ 0) , (24) 

and owing to the uniformity of the system: 

( a iL a p) = <5k,p(4. a k) , (okOp) = 5-k,p(a k a_k) • (25) 

Note that \x in Eq. (21) is the chemical potential enforcing the normalization condition 
(10). In a system without Bose-Einstein condensate, there is no need to introduce another 
Lagrange multiplier. However, as soon as the gauge symmetry is broken and a Bose-Einstein 
condensate appears, the theory acquires the new variable ^, the order parameter of the 
condensate, which satisfies the normalization condition (11). For the self-consistency of 
the theory, it is then necessary to take account of this additional normalization condition, 
which requires the second Lagrange multiplier /i . Without the latter, the theory cannot be 
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made self-consistent, and the normalization condition for the condensed particles cannot be 
guaranteed. 

It is worth emphasizing that // is the chemical potential existing for the system in both 
the gauge-symmetric and the broken-symmetry phase. At T c , the chemical potential is 
continuous, fi(T c — 0) = /x(T c + 0). One should not confuse \i with /x which is just a Lagrange 
multiplier guaranteeing the validity of the normalization condition (11), and keeping the 
theory self-consistent and the system stable. 



III. HARTREE-FOCK-BOGOLIUBOV APPROXIMATION 



We are now prepared to treat the Hamiltonian (18) with our modified HFB approxima- 
tion. For this, we introduce some notations. The momentum distribution of particles 

n k = (4a k ) (26) 

will be referred to as the normal average, contrary to the anomalous average 

o"k = (a k a_ k ) • (27) 

Summing these averages over momenta, one gets the density of uncondensed particles 

(28) 
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and the anomalous density 



Pi = 77 E ™k 



In the HFB approximation, the third-order term (22) is zero, 

=0 , 

because of (24). The fourth-order term (23), finally, becomes 

k^o v A J 

+ 77 2 9 ^k+pOpQp + ^ (o-k+paj.alp + o- k+p a_ p a p ) - - (n k+p n p + o- k+p o-*) 

k,p^0 L 

Let us define the shifted particle energies 

k 2 



(29) 



(30) 



• (31) 



£k = — + 2pg - fi , 



(32) 



where 



P = Po + Pi 



(33) 
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is the total particle density. We also introduce the notation 

A = (p + a ± )g . (34) 
Then the Hamiltonian (18) can be written in the HFB approximation as 

if H FB = ^HFB + £k4: a k + 7} J2 A (4 a -k + «-kOk) , (35) 
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k^O ^ k^O 



where 



£hfb = # (0) " ^ (2p? + ^) <?V • (36) 
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Using this Hamiltonian in the grand-canonical potential (12), we have 

dfl 



dN 

and the second derivative is given by 

d 2 Q g 



(p + Pl + a 1 )g-e , (37) 



(38) 



dN$ V ' 

Minimizing the grand-canonical potential according to conditions (14), we obtain 

£ = (P + Pi + <ri)9 , (39) 

with the stability condition g > 0. 

The Hamiltonian (35) is quadratic and can be diagonalized by means of the Bogoliubov's 
canonical transformation 

a k = u k b k + v*_ k bl k . (40) 
This brings the Hamiltonian (35) to the Bogoliubov form 

H B = E B + J2^b%, (41) 

k^O 

with the nonoperator energy term 

E B = E HFB + l 5> k -£ k ) (42) 

1 k^O 

and the Bogoliubov spectrum 

£k = v^k - A 2 • (43) 
By the Hugenholtz-Pines theorem [11-13], the spectrum must be gapless: 

lim£ k = 0, £ k >0. (44) 



Inserting Eqs. (32) and (34) into (43), we obtain for the chemical potential the equation 



V = (P + Pi ~ <ti)9 • ( 45 ) 

It is easy to check that the same chemical potential (45) follows from the Hugenholtz-Pines 
form fi = En (0,0) — Si 2 (0) by employing Green function techniques. Both real-time or 
thermal Green functions can be used, since they are just analytical continuations of each 
other [36]. Expression (45) is the standard form of the chemical potential in the HFB 
approximation (see details in the review article [4]). 

Comparing Eqs. (20), (39), and (45), we find that the chemical potential of the conden- 
sate must satisfy 

Hq = e - /I = 2<7i# . (46) 

In the broken-symmetry phase, where the anomalous average o\ ^ 0, one has /j,q ^ 0. The 
value fj, can be zero only for an ideal gas, when g — * 0, or in the Bogoliubov approximation 
[6,7], where the third- and fourth-order Hamiltonian terms (22) and (23) are neglected. But 
in general, the Lagrange multiplier /io is non-zero, thus making the theory self-consistent. 
Using the chemical potential (45) in Eq. (32), we have 

& = £ + (47, 
and the spectrum (43) takes the form 



with the sound velocity 



r k)2+ G9 ■ (48) 



c=J-. (49) 



In the long-wave limit k — > 0, the quasiparticle energy behaves like ck, thus being gapless, 
as it should. 

For the normal average (26), we find 

n k = (ul + vl) / k b + vl , (50) 

while for the anomalous average (27), we get 

<7 k = (1 + 2f£)u k v k , (51) 

where the momentum distribution of bosonic quasiparticle excitations 

4 b - (bW) = (52) 

can be written in the form 
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/ k b = ^(coth^-l). (53) 



2T 

The coefficient functions of the Bogoliubov transformation (40) satisfy 

K - v l = 1 » M k^k = - , (54) 



2e 



k 



and 



_ . del + A 2 £ k 
4 + ^ = — = (55) 



^k ^k 



and read explicitly 



M k = ^- = —7T — , (56) 



2e k 2£l 



and 



t4 = J g^ zf . = ^_ 



2e k 2e k 
In this way, for the normal average (26), we obtain 



n k = 2^ coth— - - , (58) 

while the anomalous average (27) becomes 

(j k = - — coth— . (59) 

The anomalous average (59) is important as compared to the normal average (58). For 
this purpose, consider these averages as functions of e k in the range 

< e k . (60) 

At low momenta or energies, such that £ k <C A and e k <C T, the asymptotic behavior of the 
normal average (58) is 

TA ATI 
nk "^T + l2T + 2A " 2' (61) 



k 

whereas the anomalous average (59) is 

TA A elA 
ak ~~ll ~ 12T + 720T3 • ( 62 ) 

From here it is evident that in the long-wave limit the anomalous average is of the same 
order of magnitude as the normal one, n k ~ |cr k |, only their signs are opposite. 
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In the short-wave limit, when e k 3> A and £k » I 1 , the asymptotic behavior of the 
normal average is 

"k^fAV-fAV + e-^, (63) 



,2eJ \2e k 
whereas that of the anomalous average is 

<7k ~ - ^- (l + 2e~^) . (64) 

In this limit, the magnitude of the anomalous average becomes much larger than the normal 
one. 

Summarizing, we conclude that, in the above two limits, we have 

1 0k | —n k (k — > 0) , |<7k| > n k (A; — > oo) . (65) 

Consequently, the anomalous average is always important. 

In the large- volume limit, the sums (28) and (29) can be calculated as momentum inte- 
grals, leading to the density of uncondensed particles (28): 



Pi = n 



d 3 k 



2 J (2tt) ; 



£k 



(66) 



while the anomalous density (29) becomes 

d 3 k A /fe k 



(2vr) 3 2e k 

Taking into account relation (34) and using the notation 



^-y (2^^ coth T- (67) 



we obtain 



/• d 3 k g fie^ . , 

a = Jj2^2T k C ° th — ' (68) 



A = -M-, (69) 
1 + a 

so that the anomalous density (67) can be represented in the form 

ai = ~ TT^ ' (70) 

The quantity (68) is ultraviolet-divergent. This divergence is caused by the modeling 
of the short-range interaction in the dilute-gas approximation by a 5-function potential 
(1). There are two ways of removing this divergence. One may either use a more realistic 
interaction potential V(x), whose Fourier transform goes to zero for k — > oo, for instance 
Gaussian-type potentials [37-39]. More efficiently, with the same physical consequences in 
the dilute limit, one renormalizes the 5-function potential, replacing it by the scattering 
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matrix obtained from the Lippmann-Schwinger equation. This simply renormalizes the 
coupling constant g to the renormalized defined by 



9R g J (27r) 3 kV2m-iO' 

Indeed, the relation to the scattering length in Eq. (2) is really only valid for the renormal- 
ized coupling constant gn, which is usually not mentioned for brevity, until the ultraviolet 
divergences appear. This procedure is standard for eliminating divergences in calculating 
the ground-state energy [2-4] as well as the anomalous averages [40-42]. Keeping in mind 
such a correction, we may consider the quantity (68) as finite. 

It is easy to check that condensation proceeds in a second-order phase transition. When 
p — > 0, p\ — > p, we see that A — > and — > k 2 /2m. From Eq. (66) we find the 
condensation temperature 



A c 

m 



P 



am 



2/3 

(72) 



which coincides with that of the ideal gas, as it should be for a dilute gas in mean-field 
approximation. 

IV. ZERO TEMPERATURE 

The above equations can be calculated explicitly in the zero-temperature limit, where 
the density (66) becomes 



Pi 



This can also be represented as 



1 r d 3 k 



2 J (2vr) 3 V^k 
Substituting Eqs. (47) and (48), we find 



3tt 2 

For the fraction of condensed particles 



P. = ^ ■ (75) 



P P 

we obtain 

■ ■• 

3tt 2 



n = ^ = l-n 1 = l-^, (76) 



(mca) / n \ 

no = 1 - ^"o-T" (P a = X ) ' ( 7? ) 
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where a is the mean interparticle distance. 

The integral (68) becomes, after the renormalization of g according to (71), 



a 



2 J Ve k k 2 J (2tt) 3 ' 
where is sloppily replaced by g, as usual in such calculations. The value of a is 

1 



a = ^ m 2 cg , (79) 

so that the anomalous average (70) becomes 

p m 2 c g 

o-i = — — . (80) 

7r z — m A c g 

This way of removing the ultraviolet divergences is quantitatively exact for calculating quan- 
tities depending only on the scattering length a s in one-loop approximations [4]. For strong 
interactions, the shape of the interaction potential becomes important [4,37-39]. In what 
follows, we shall consider the anomalous average (80) for arbitrary g, but keep in mind that 
at large g this expression for o\ can be only qualitatively correct. Being proportional to 
p , the anomalous density a\ tends to zero for p — > 0. However, in a dilute gas at low 
temperature, p is close to p as follows from Eq. (77). 

In further calculations, it is convenient to introduce the diluteness parameter 



5 = p 1/3 a s = ^ (81) 
a 

and to introduce the reduced dimensionless sound velocity 

c = cam . (82) 
With this notation, the fraction of condensed particles is given by 



c 3 



no = 1 - — 2 , (83) 



and the fraction of uncondensed particles by 



?,3 



ni = - 2 . (84) 



Combining Eq. (79) with relation (2) for the renormalized g, we obtain 



4 r 

a = - - c5 . (85) 

7T 



The anomalous average (80) is then 



<?i = ~~pr7 no . (86) 
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To define the above quantities (83)-(86) as functions of the diluteness parameter (81), 
we have to know the dependence on 5 of the reduced sound velocity (82). For this purpose, 
we use Eq. (69) in the form 



47r 2 £ / n 



(87) 



7r — Ac5 \ma 2 , 

and recall that A = mc 2 according to Eq. (49). Thus Eq. (87) becomes 

85c 3 - 3ttc 2 + 12vr 2 5 = , (88) 

which fixes the reduced sound velocity as a function c(S) of the diluteness parameter. There 
exists a positive solution of Eq. (88) for c in the (^-interval < 5 < 5 C , limited by the critical 
value 

1 /TlV/S 

5 C = - I - = 0.253873 . (89) 



4 \3J 

Remarkably, there exists an interesting relation between the anomalous average (86) and 
the density of uncondensed particles p±. We take the identity following from A = mc 2 : 

Ama 2 = (mac) 2 = c 2 , (90) 

and substitute on the left-hand side A from Eq. (87) to obtain 

c 2 = ^ . (91) 

7T - 4cS 

Inverting this with respect to 5 yields 



nc 2 



(92) 



4(7r 2 n + c 3 ) 
Substituting Eq. (92) into Eq. (86), we obtain 

a l = 3 Pl . (93) 

That is, the anomalous average is three times larger than the normal one. This emphasizes 
once more that at low temperatures the anomalous averages are important. 

At asymptotically small diluteness parameter 5 <^ 1, the reduced sound velocity (82) in 
Eq. (88) behaves as 

c ~ 2^ S 1 / 2 + l 4 S 2 + -^L 5 7 ' 2 (5 « 1) , (94) 

3 9^/7T 

while the condensate fraction (83) has the expansion 

Keeping only the first term in the expansion (94), we obtain the Bogoliubov sound velocity 
cb = ^4:7ipa s /m. Retaining on the right-hand side of Eq. (95) the first two terms, we obtain 
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the Bogoliubov depletion formula [6-9]. In the Sec. V, we show that the ground-state energy 
in the limit 5^0 also gives the known Bogoliubov expression. Thus, in the limit of the small 
diluteness parameter 5, our equations have the correct asymptotic behavior, reproducing the 
results of the Bogoliubov approximation. 

When the diluteness parameter approaches the critical value (89), the sound velocity 
c = c/am has the expansion 

1 /3 

c~c c -2^(5 c -5) 1 / 2 + 8(|) (5 c -5) (5 - 5 C - 0) , (96) 

with 

c c = (3tt 2 ) 1/3 = 3.093668 . (97) 
And the condensate fraction is given by 

n ~ 6 Q 1/6 (5 C - 5) 1 ' 2 - 20 (^)'* (S c -5) (S^5 c -0). (98) 
At the critical depletion (89), we have 

c=^^, n = (5 = 5 C ). (99) 
am 

The disappearance of the condensate fraction for 8 > 5 C at zero temperature is a signal 
for a quantum phase transition. Here this transition occurs as a function of the diluteness 
parameter 5. To display the critical behavior of physical properties, it is convenient to 
introduce the relative distance variable from the quantum critical point 

r=-^. (100) 

Oc 

Then the reduced sound velocity behaves in the vicinity of the critical point as 

^— ^ * -r 1/2 + 2 -r (r « 1) . (101) 

c c 3 

For the condensate fraction we obtain 

n ~ 3r 1/2 - 5r (r < 1) , (102) 

implying that the critical exponent (3 for the order parameter \1/ is 1/4. 

The overall behavior of the reduced sound velocity c = c(5) and of the condensate 
fraction no = Uq(5) as functions of the diluteness parameter 8 are shown in Figs. 1 and 2, 
respectively. The function c(S) in Fig. 1 is calculated from Eq. (88). Substituting this c(S) 
into Eq. (83), we get the condensate fraction n (5) plotted in Fig. 2. 
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V. THERMODYNAMIC CONSISTENCY 



As a final important point we convince ourselves that our self-consistent bosonic HFB 
approximation involving the new Lagrange multiplier /i is consistent with the thermody- 
namic formalism, in contrast to previous attempts to generalize the HFB approximation to 
condensed Bose systems. For the discussion of the inconsistencies we mention once more the 
work of Hohenberg and Martin [21] and the subsequent papers [17-19,30-32,43]. The most 
recent and very clear analysis of the problem was given in the review article by Andersen 
[4]- 

The grand-canonical potential in Eq. (12) and the associated Hamiltonian H in Eq. 
(13) make it possible to find all thermodynamic properties of the system. The free energy 
(9) is connected with the grand potential (12) through the relation 

Q — F — /i A - fxN . (103) 

The Lagrange multipliers uo and /i should not be confused with the standard chemical 
potential defined for the system without condensate. In the presence of a condensate, its 
role is played by what we may call effective system chemical potential, denoting it by Ji. In 
the condensed phase, it is given by 

p, = li n Q + ii, (104) 

where no = N /N. With this definition, the grand-canonical potential (103) satisfies the 
usual thermodynamic relation 

Q = F-fiN. (105) 

Its differential satisfies 

oin = -SolT - PalV - Ndfi , (106) 
where S is entropy and P pressure. For the free energy F, this implies 

dF = -SalT - PalV + fialN . (107) 
Thus the system chemical potential fi is given by the derivative 

From these expressions, we may calculate all thermodynamic properties of the condensed 
Bose gas. Thus our self-consistent HFB approximation obeys the standard thermodynamic 
formalism. 

As an illustration, let us verify the consistency of the above two definitions (104) and 
(108) for the system chemical potential. To be explicit, consider a dilute Bose gas at zero 
temperature, representing the universal terms [4], which are asymptotically exact for all 
short-range interaction potentials with scattering length a s . 

At zero temperature, the free energy coincides with the internal ground-state energy, 
F — E. Then Eq. (108) reduces to 
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A = h^7 (T = 0) . (109) 



The internal energy is, by definition, 

E = (H) + fiN , (110) 

where if is the grand-canonical Hamiltonian (13). The average (H) is given in the HFB 
approximation by 

(H) = E B = E HFB + Vj^(e k -^), (111) 
where, according to Eqs. (42) and (19), 

£hfb = # (0) - ^ (2p? + <r?) = iV (I Po „ - „ - - ^ (2p? + <r?) . (112) 

For the explicit calculation of the integral in (111) in the dilute gas, we perform the usual 
subtraction implied by the renormalization of g in Eq. (71). The subtracted integral becomes 



d 3 k 



1 r d 6 i 

2 J CM 



£k - £k + 



3 4 

m c 



(2tt) 3 [ k 



C m 4 5 



Employing the results of Sec. IV, we find for the internal energy (110), which at T = is 
the ground-state energy, 

e^JHU + J^M. (1 i4) 

2V V 15^ J v ; 

This is the expression derived by Lee et al. [44,45] for a hard-sphere potential. It is uni- 
versal in the sense that it applies to any short-range potential with scattering length a s [4]. 
Differentiating (114) with respect to N we obtain, according to (109) and using the fact that 
F = E at T = 0, the effective chemical potential 

~^pg(l + ^=5^ , (115) 

where we have inserted the derivative 85/ dN = 5/3N. The same result is obtained from 
(104), showing the self-consistency of our HFB approximation. 

It is important to emphasize that the self-consistency in our approach has been achieved 
by accurately taking into account all normalization conditions, which required the intro- 
duction of the additional Lagrange multiplier /x . It is due to the latter that our HFB 
approximation is completely self-consistent and displays no gap in the spectrum. Without 
Ho, we would plunge back to the known problem of the standard HFB approximation, which 
is not self-consistent and possesses an unphysical gap in the spectrum [17-19,30-32,43]. 

In conclusion, we have presented a new bosonic self-consistent Hartree-Fock-Bogoliubov 
approximation. It is derived from a variational principle and preserves all symmetries of 
the Hamiltonian. At the same time, it is gapless in the condensed phase, thus solving an 
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old outstanding problem of Bose systems. We did not invoked any unjustified tricks as in 
previous attempts with the same goal, such as omitting anomalous averages, and avoided 
any mismatch of approximations by adding additional phenomenological terms to remove the 
gap. This became possible by accurately taking into account two normalization conditions. 
Thus our HFB approximation is conserving, gapless, and self-consistent. 
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FIGURES 




FIG. 1. The reduced sound velocity c = c(<5) as a function of the diluteness parameter 5, which 
is given by the solution of Eq. (88). 
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FIG. 2. The condensate fraction no = no(5) as a function of 5, obtained from Eq. (83). 
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Figure Captions 



Fig. 1. The reduced sound velocity c = 6(6) as a function of the diluteness parameter 6, 
which is given by the solution of Eq. (88). 



Fig. 2. The condensate fraction n = n (6) as a function of 6, obtained from Eq. (83). 
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